#ifndef SOURCESTENCIL_H
#define SOURCESTENCIL_H

#include "Core/Tensor.h"
#include "Embedded/EmbeddedGrid.h"

template <int Dim>
struct SourceStencilParam
{
  enum class OperatorType { Elliptic, Convection, Others };
  OperatorType otype;
  // For convection type
  const Tensor<Real, Dim+1> *velField;
};

template <int Dim>
class SourceStencil
{
public:
  using rVec = Vec<Real,Dim>;
  using iVec = Vec<int, Dim>;
  using EmGrid = EmbeddedGrid<Dim>;
  template <class T> using vector = std::vector<T>;

  template <class T_DiffOp>
  SourceStencil(const EmGrid &aGrid, const T_DiffOp &diffOp, const SourceStencilParam<Dim> &param);

  void dump(std::ostream &os) const;

protected:
  // Place the center of the box in the node i s.t. norm(anchor - i, 0) \in [m, m+1).
  template <int n>
  bool generateLatticeOnRing(const iVec &p, int m, const iVec &q0, std::vector<iVec> &outLoc) const;

public:
  struct Stencil {
    iVec         destIdx;
    vector<iVec> locations;
    vector<int>  isBdry;
    vector<Real> weights;
  };
  vector<Stencil>  vStencil;
  const EmGrid     &emGrid;
};

#endif //SOURCESTENCIL_H
